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1 Introduction 

It is well known that computing the eigenvalues of a self-adjoint bounded or differ- 
ential operator H acting in an infinite dimensional space Ti is not straightforward. 
The standard truncation (or Ritz) method is to take a large finite-dimensional test 
function space C with an orthonormal basis {0j,}r=i ^^cl to compute the eigenvalues 
of the n X n matrix A given by 

Ar,s = {H(f)r,4>s)- 

Unfortunately one often encounters the problem of spectral pollution: the ma- 
trix may possess eigenvalues which are unrelated to any spectral properties of the 
original operator. This is an important problem in several areas of applied math- 
ematics, as is illustrated by recent papers, [H |2l El IH ISl UHl UH UIl 120] • There 
are methods of detecting and resolving such problems, but they are not easy to 
understand and for that reason are not well appreciated. Our main result in this 
paper is the proof that, at least in some situations, the geometrically motivated 
method of Davies, |7j is mathematically equivalent to the method of Zimmermann 
and Mertins, [20] which is based on Lehmann's method ^3] and its extension by 
Goerisch PP. We also provide various other ideas and examples to motivate the 
method and show how well it works. 



2 Formulation of the Problem 

We describe the context of the paper at a rather abstract level. The assumptions 
are designed to apply to a range of problems, including those involving second 
order elliptic operators, as will become clear in our applications. 

Let A be a bounded, self-adjoint, hnear operator acting on a Hilbert space H, and 
let P be the orthogonal projection onto a closed linear subspace £ of TC. The 
truncation of A to £, denoted Ac, is defined to be the operator PAP restricted to 
C (Similar definitions may be made for unbounded A and for quadratic forms.) 
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If £ is sufficiently large, one might hope that the spectrum of Ac is close to that 
of A. The Rayleigh-Ritz, or variational, theorems establish precisely this for the 
part of the spectrum of A outside where a (resp. /5) is the smallest (resp. 

biggest) point in the essential spectrum of A. 

Within [a, f3] the situation is far less pleasant. One says that A is a point of spectral 
pollution of A for the increasing sequence of subspaces with union dense in H 
if there exist A„ G Spec{Ac„) such that A„ ^ A as n ^ oo but A ^ Spec(y4). We 
will see that spectral pollution is commonplace, but the following two elementary 
lemmas sometimes enable one to determine Spec (A) from its truncations. 

Lemma 1 // ||(/ - P)AP\\ < e and A G Spec(A£) then 

(A - e, A + e) n Spec(A) ^0. 

Proof One considers A as a perturbation of PAP + {I — P)A{I — P). 
Lemma 2 // ||/|| = 1 and \\Af - \f\\ < e then 

(A-£,A + £) n Spec(A) ^ 0. 

Proof This uses the resolvent norm estimate 

\\{A - Xiy^W = dist{A, Spec(A)}-^ 

Lemma 2 is originally due to D. Weinstein. Methods of estimating Spec (A) starting 
from a subspace C may be divided into two general categories. Prior methods 
involve proving that if is a particular increasing sequence of subspaces then the 
spectrum of Ac„ converges to that of A. There should also be an explicit error 
bound, or at least some information about the rate of convergence. When applied 
to differential operators such methods often involve finite element subspaces and 
Sobolev space embedding theorems. They have two disadvantages. The ffist is 
that the convergence proofs are often very difficult and assume conditions which 
are specific to the particular applications. The other is that the error bounds are 
often far larger than the true errors. Nevertheless valuable and deep theoretical 
results can sometimes be obtained using such methods, [121 US] • 

Posterior methods, on the other hand, involve choosing a subspace C and carrying 
out computations which provide rigorous bounds on the location of Spec(yl). If 
these are good enough one stops, but if they are not one chooses a different or 
larger subspace C The problem here is that there is no prior guarantee that useful 
bounds will ever be obtained. Nevertheless in practice the method is often a very 
good one. Besides the result by D. Weinstein (Lemma 2), the more refined spectral 
bounds by T. Kato [llj, Aronszajn-A. Weinstein ^^1, Lehmann jHj and Goerisch 
fall into this category. Two particular posterior methods have recently been the 
focus of some attention. One uses the second order spectrum, introduced in [71ITH]. 
and is the subject of the companion paper by Levitin and Shargorodsky, jTHj. This 
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paper studies the other, which depends upon a systematic study of how certain 
residuals depend upon a real parameter, 7j. We describe a close relationship 
between this method and that of Zimmermann and Mertins, [20] . 

We next set the technical context of the paper. Let if be a non-negative unbounded 
self-adjoint operator acting on 7i, and suppose that ^ Spec{H). Now suppose 
that K is a second non- negative self-adjoint operator, possibly with a different 
operator domain, but satisfying a quadratic form relative bound. In other words 
if Qh and Qk are their two quadratic forms then Dom^Qn) C Dom((5x) and 

< QkU) < aQnif) 

for all / e Dom{QH) = Dom(if^/^). (In fact Qx need not be closed, and hence 
need not be determined by an operator K, in the analysis below.) 

The eigenvalue problem XHf = Kf may be recast in the weak, or variational, 
form. In other words one seeks /, A such that 

>^QH{f,g) = QK{f,g) (l) 

for all g G Dom(Q/^^). 

In many cases the operator K is of the following form. There exists a closed 
operator D : Ti ^ M. with operator domain Dom(Qj^) such that H = D*D and 
K = D*MD, where M is a non-negative bounded operator acting on Ai. Under 
our hypotheses Ker(D) = and Ran(D) is closed. There is no loss of generality 
in assuming that Ran(D) = Ai, hj simply reducing the size of Ai and restricting 
M. In this case the issue is to find / such that 

X{Df,Dg) = {MDf,Dg) 

for all g, or equivalently to find the eigenvalues of the bounded operator self-adjoint 
M. 

Many non-negative differential operators H of order 2n can be written in the form 
H = D*D, where D is a differential operator of order n. Indeed, in a formal 
sense the assumption is always satisfied with D = H^^^ and M = H~^^^KH~^/^ . 
This observation is rarely useful, because H^^/'^ are usually not easy to describe 
in explicit terms. If H and K are both differential operators, and K is of lower 
order than if, then one often finds that M is compact. Its essential spectrum is 
then simply {0}, and the pollution problem does not arise. 

Before continuing, we mention that for Sturm-Liouville equations one can often by- 
pass the problems in this paper by solving the initial value problem and computing 
the Prufer angles of the solution. This can enable one to locate eigenvalues in a 
gap of the essential spectrum. The method is not available in higher dimensions. 

3 Reality of the Phenomenon 

Experience shows that spectral pollution is not just a theoretical problem. It 
occurs on a regular basis for operators of physical importance, such as truncations 
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of periodic Schrodinger operators to finite regions; see for example [Hj. In this 
section we demonstrate how bad the phenomenon can be. Later in the paper we 
produce more reahstic examples. 

The following well-known argument (see e.g. JHl) shows that spectral pollution 
can even occur in an extreme manner for basis functions which are finite linear 
combinations of exact eigenfunctions, and which are therefore in the operator do- 
main. 

Lemma 3 Let (A^°^)„gN and (A^"^^)„gN be two (possibly constant) sequences of 
eigenvalues of problem (QJj, converging to limits X^^^ and X^^\ respectively, sat- 
isfying A^^-* < A*^"^-*. Moreover, let jj, satisfy A'-'^-* < jj, < A^"'^^ Then, there exists a 
sequence {ipn) of trial functions which all are finite linear combinations of exact 
eigenfunctions of ([TJj, such that 

QH{i'n,i'm) = Snm, QxO'ny'ipm) = f^^nm, (2) 

i.e., /i occurs as the only eigenvalue of the Ritz method, for every dimension. 

Proof Without loss of generality, let < /i < A^^^ for all n E N, and let (/^°^) 
and (f^^) denote sequences of eigenfunctions of (jT} corresponding to (A^°^) and 
(A^^^), respectively, such that 

^H\Jn 1 Jm J '^nmj mK\Jn i Jm J '-'nm- 

Finally, choose G [0, 27r) such that 

(cos^a„)A^°^ + (sin^a„)A^-^^ = /i, 

and 

ipn ■■= (cosa„)/^°) + (sina„)/^^^ 
(j2I) is then quite obvious. 

In the same way, also arbitrary sequences (/in) between (X^^) and (A^^^) can be 
shown to occur as spurious eigenvalues. 

If we make the weaker assumption that A*-'^-' and A*-^-* are in the essential spectrum 
of M, then a similar argument can be followed. One can ensure that Dip^ lie in 
the domain of all powers of the operator M. 

One may criticize the above example on the grounds that the sequence {ipn) of basis 
functions constructed above may be far from complete. This gap can be closed as 
follows, at least if problem has a Qj^-orthonormal and complete system of 
eigenfunctions: 

In this case, let (/^^^)n6N denote a third sequence of eigenfunctions such that 
ifn^)^ ifn^)^ ifn^) together are Q/^-orthonormal and complete. Then, with 
(q;„) as above, define 

V?! := (cosai)/f^ + (sinai)/^^\ 
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and for n > 1: 

■= {- sin an) fi°^ + {cos an) fi^\ 
fAn-i ■= (cos aan) /in + (sin a2n) /in , 

!fAn ■■= (cOSa2n+l)/in+l + (sinasn+O/in+l, 

Then, (v5n)nGN is (5j/-orthonormal and complete, and the Ritz method (with basis 
functions of dimension provides yU as a ~ ^-fold eigenvalue. 



4 Methods of Avoiding Pollution 

There are various ways of identifying and avoiding spurious eigenvalues. Sup- 
pose that (Dom{QH),QH) is a Hilbert space, with continuous embedding into H. 
Let (£„)„gN denote an increasing sequence of subspaces, with UnGN-^n dense in 
Dom((5//). Let {fin, Xn)neN be a sequence of approximate eigenpairs to the problem 
(|H) obtained by the Ritz method: 

Xn e £„, HnQH{Xn, q) = QK{Xn, q) for all g G £„. (3) 



Lemma 4 Suppose that /i„ — >■ G R (or that a subsequence converges to n), with 
/i not an eigenvalue of (OP; is therefore either in the essential spectrum or a point 
of spectral pollution. 

Then, if we normalize Xn by 



\Xn\H ■= \/QH{Xn,Xn) =1, (4) 

we have 

Xn ^ weakly in {Dom{Q h) , Q h) ■ (5) 

Note If in particular the embedding Dom((5/f) 7^ is compact, we therefore 
obtain x„ ^ in {H, ||.||), i.e., 

llXnII , Q_ 



Since in practice one often performs the Ritz method for various different dimen- 
sions n, one has the chance of getting evidence for (jHl) by evaluating the quotient 
for various n. 

Proof We show that each subsequence of {xn) has a subsequence converging weakly 
to zero. Indeed, (jl} implies that a subsequence (Xn*,) of a given subsequence 
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converges weakly in Dom^Qn) to some x ^ Dom((5H). Since Qk > is bounded 
in (Dora{QH),QH), this implies QxiXrik^g) ^ Qk{x^9) for each g e Dom(QH)- 
Consequently, 

f^nkQuixnu^g) - QK{xnt,,g) ^^Qh{x, g) - Qk{x, g) (7) 

for each g G Dom(Q//). By Q, the left-hand side of ((Tj) is zero for g E Cn and 
nfc > n. Therefore, ((Tj) implies 

fj-Qnix^g) = QKix,g) for all e |J 

nGN 

and thus, by density, for all g G Dom^Qn)- Since /i is no eigenvalue of ((H), the 
assertion x = follows. 

One may find the eigenvalues of a bounded operator M in a gap of its 

essential spectrum by considering instead the operator = (M — pl)~^, for some 
p G (a, The eigenvalues of M in (a, /5) are associated to eigenvalues of which 
lie outside its essential spectrum; the latter may therefore be found by using the 
Rayleigh-Ritz procedure. An apparent disadvantage of this approach is the need 
to find an expression for the operator A^, but this may be avoided, (20] . 

In this paper we explore a geometrically motivated method of finding eigenvalues, 
which traces its origins back to Kato, ^1], and Weidmann, [HI Cor. 6.20]. One 
of our main goals is to describe its close relationship to ideas of Goerisch and 
Lehmann, and to developments by various other authors, [H UHl ED! • Suppose that 
the self-adjoint operator M acting on the Hilbert space Ai has a single eigenvalue 
/i which lies within a gap in its spectrum. If 

|A — yu| < min{|A — a\, \X — 

then one may apply the Rayleigh-Ritz method to the operator B = (M — A/)^ for 
any A close enough to fi. The smallest eigenvalue of B will then be (A — /i)^, from 
which fi can be calculated (up to a square root ambiguity). The idea is to choose 
the value of A to maximize the numerical efficiency of the method. One's intuition 
is to take A as close as possible to the eigenvalue p, but we shall see that this is 
not the correct strategy. 

As in iTj, we investigate the A-dependence of the procedure above by introducing 
the function 

F(A) = min{||M/ - A/|| : / G and ||/|| = 1}. 

If one could evaluate F numerically then the formula 

F(A) = dist{A,Spec(M)} (8) 

would enable one to locate the spectrum of M. To find the eigenvalues of M one 
needs to combine this with some theoretical input which distinguishes between the 
point spectrum and essential spectrum of M. 
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Following [7j, we introduce the approximate functions 

F„(A) = min{||M/ - A/|| : / G and ||/|| = 1}. 

These may be used to obtain useful spectral information, as follows. 

Lemma 5 If M is bounded then for any e > there exists such that n > 
implies 

Fis)<Fn{s)<F{s)+e (9) 

for all s G R. 

Proof It follows from their definition that Fn{s) decrease monotonically to F{s) 
as n — > oo for each s G R. It also follows from its definition that 

|F„(A) < |A-/i| 

for all A,/i G R. Equivalently |F'(A)| < 1 for almost all A G R (in the sense of 
Lebesgue). Since the functions F„ are equicontinuous they must converge uniformly 
on compact subsets of R. Suppose that n > implies © for — ||M|| < s < \\M\\. 
If s > ||M|| then 

F{s) < Fn{s) 

< F„,(||M||) + .-||M|| 

< F(||M||)+£ + .-||M|| 
= F{s)+e. 

A similar argument applies if s < — ||M||. 

Lemma 6 Given (0), let A be an eigenvalue of M and let there be no other point 
of Spec(M) in [A — Ae,X + Ae]. Then Fn has a local minimum at some point 
a E [\ — e, X + e]. If n = a — Fn{a) and u = a + Fn{a) then 

X-2e<fi<X<u<X + 2e. 

Proof Since F{X) = we see that < Fn{X) < e. However if e < |s — A| < 2£: we 
have 

e < |s - A| = F{s) < F„(s). 

Hence the minimum of F„ within the interval [A — 2e,X + 2e] must occur for a 
within [X — 6, X + e]. Moreover < -F„(cr) < e. The second statement of the lemma 
follows immediately. 

The functions Fn{-) may be computed as follows. One introduces 

F^,„(A) = min{||P„(M/ - A/)|| : / G £„ and ||/|| = 1} 

where m > n. U m = n then this achieves no more than finding the eigenvalues of 
the truncation of M to £„ (i.e. applying the Ritz method), but if m ^ it may 
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be far superior. One may also compute F„(A) exactly by replacing above by 
the finite-dimensional subspace 

The projection in the argument below is then replaced by the projection P„ 
onto £„. Whether or not this Lanczos process is easily implemented depends upon 
the example. 

The evaluation of Fm,n{^) is straightforward. One has the formula 

where G(m, n, A) is the smallest eigenvalue of the n x n matrix 
B{m, n, A) = PnMP^MP^ - 2AP„MP„ + A^P^. 

In practice one finds potential eigenvalues A of M by using the Ritz method. One 
then examines Pm,n(A) for those values. If it is sufficiently small then one has good 
evidence that there is indeed an eigenvalue near the candidate value, while if it is 
not small then one has no such evidence. Using only a finite number of matrix 
entries of the operator M no more could ever be established, but this method is 
far better than simply using the Ritz method without any checks. 

The above provides evidence for the existence of an eigenvalue near a specified 
number. One can also obtain evidence for its multiplicity, or more precisely for 
the number of eigenvalues in a very short interval, by an elaboration of the above 
ideas which is described in |7]. 



5 Bounding the eigenvalues 

In [7j we discussed how to obtain rigorous bounds on the spectrum of M (or of 
problem ([T]), respectively) from the functions Fn{-), thus avoiding the problem of 
spectral pollution. In this paper we develop the ideas further. If the subspace 
is large enough then Fn{X) has local minima near the eigenvalues of H, and its 
derivatives are close to ±1 on either side of the eigenvalues. The following exact 
formula for the derivative may be useful, for example if applying Newton's method. 
We have 

= 

where Gn(A) is the smallest eigenvalue of the n x n matrix 

P„(A) = PnM*MPn-2XPnMPn + X^Pn (10) 

= Dn-2\Mn + \^I (11) 
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where the second hne is understood as acting in Ln and M„ = -/^£„- If /n(A) is the 
normahzed eigenvector corresponding to the eigenvalue G'„(A) then it is well-known 
that 



A computation of the graph of the function F^i^ is bound to be slow if n is large. It 
is best therefore to attempt this only for reasonably small n, in order to understand 
approximately where the spectrum is located. One then passes to larger values of 
n for which F„(A) is only computed for carefully selected values of A. 

We next discuss the computation of the graph of -F„(') for reasonably small n. The 
procedure below is guided by the fact that it is optimal if Dn = M^, or cquivalcntly 
if (/ - Pn)MPn = 0. It should therefore be close to optimal if Dn - = 
PnM{I — Pn)MPn is Small. The starting point is to find all of the eigenvalues A^ 
of M„ and then to compute G(Aj.) for each r. If G(Aj.) is not sufficiently small for 
some r then A^ is rejected as spurious. Let the remaining eigenvalues be listed in 
increasing order. The corresponding normalized eigenvectors gr- of are used as 
starting vectors for an inverse power iteration to find = f{\)- As a check it 
should be verified that fr and Qr are approximately equal, and that the set fr is 
approximately orthonormal as r varies. 

Now assume that one wants a list of values of Fn{-) over a set of points covering 
the interval (A^, A^+i), in order to plot a graph. This is done separately on either 
side of the midpoint 7^ of the interval, and we consider the left hand half. One 
chooses some subdivision A^ = di < o"2 < ••• < o"m < 7r and then computes G{(Jr) 
in order. At each stage the starting vector for the inverse iteration should be the 
eigenvector for the previous stage. In the right hand half interval one deals with 
the points in decreasing order. 

The reason for using different procedures on either side of the midpoint becomes 
clear if one considers the exact case, in which = M^; the eigenvector is then 
constant in each subinterval and changes discontinuously as A increases through 
7j.. In the inexact case one should expect /(A) to change rapidly near A = 7,., and 
it is therefore advisable to start the iteration with the vector /(A') + /(A"), where 
A' and A" are values already computed on either side of 7^. 

One also needs to evaluate F„(A) for large values of n at individual points A satis- 
fying Xr < X < Xr+i- Once again this may be done using inverse power iteration, 
but using fr + fr+i as the initial vector. 

We finally discuss how the computed values of F„(A) are used to locate the spec- 
trum of M. We use 



^^(A) 



K(A)/„(A),/„(A)) 
2A-2(M„/„(A),/„(A)). 



(12) 
(13) 



Spec(M) n [A 



F„(A),A + F„(A)]^0 



for all A e R, which follows from the spectral theorem, and the obvious 



Fn{X) > F(A) 
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for all A G R. Unless F„(A) and F{\) are approximately equal, the above provide 
little information. One has two choices: 

(i) One can supplement the computation of -F„(A) by prior information about the 
location of the spectrum of M. 

(ii) One can examine the graph of Fn{X) to see whether its form is consistent with 
the assumption that it is approximately equal to F{\), using (jS)). If not then one 
can increase the value of n, or revert to (i). 

Of course the former procedure is to be preferred, but the latter is likely to be 
easier. We illustrate both in the example discussed in the next section. 

For the rest of the paper we make the following hypothesis. 

(S) We have > F{s) and < 1 for all s G R. 

This can be ensured by replacing F„ by 

F„(s) = F„(s)+£e-^"(^) 

for any e > 0. Increasing Fn{s) slightly makes the spectral bounds which we prove 
slightly weaker, but does not involve an essential loss of generality. We will use the 
following consequence of the hypothesis (S). 

Spec(M) n (A - F„(A), A + F,(A)) ^ 0. 

Our goal is to determine those eigenvalues of M which lie in a given interval {a, (3); 
a and /3 might sometimes be known points of the essential spectrum of M. Let 
(Jr, 1 < r < R he an increasing sequence of numbers lying in the interval (a,/?); 
they are our first approximations to the eigenvalues of M within the interval. They 
might be the non-spurious eigenvalues A^ of M„ already described, or they might 
be the local minima of the function -F„(A), obtained by solving F^{\) = using 
the bisection method. In the simplest cases these two methods of choosing cxr yield 
approximately the same values. If this fails one can still proceed, using the second 
choice of ar just mentioned. Putting fir = cr^ — Fn{ar) and Ur = (Tr + Fn{(Jr) for 
1 < r < R, we have 

Spec(M) n {fir, 7^ 

for such r. We also put = a and /i_R+i = p. Further progress depends upon the 
following hypothesis. Following (i) above, this might have been proved, possibly 
using a homotopy method. Or following (ii), it might simply be motivated by the 
numerical evidence, and tested by the consistency of the consequences drawn using 
it. 

(H) We have 

Vq< fll<Vl< fl2< ■■■ < fJ'B. < < l^R+l- 

The intervals {ur-i, fir) do not meet Spec(M) for 1 < r < R+ 1, while [fir, Vr] each 
contain a single eigenvalue nir of M for 1 < r < R. 
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Situations in which the operator M has eigenvalues which are too close to be 
resolved using the subspace can often be handled using a modified procedure, 
or more simply by increasing n. 

We describe an iterative procedure for obtaining steadily more accurate enclosures 
of the eigenvalues. The following lemma provides the starting point for the it- 
eration. The constructions are most easily understood by plotting the graph of 
a typical function F„. See Figure 1 for an example in which a — 0, (3 — 1 are 
eigenvalues of infinite multiplicity, and R—1. 

Lemma 7 Assuming (S) and (H), we have 

for 1 < r < R, where v'^ and ii'^ are defined by solving the following equations. 

-^n(^r) — ^r— 1 

Z/' = Sr + Fn{Sr) 

Fnitr) = fJ'r+1 — tr 
l-^r Fnitr) 

Proof Let 1 < r < R. The function 

= Fnis) - S + Vr-i 

satisfies ^'(s) = F^{s) — 1 < for all s, 0(^'r-i) = > and 

0((Tr) = Fniar) — (7^ + I'r-l = — < 0. 

Hence Ur-i < Sr < Cr- We next observe that 

Fni^r^ Fni^r^ ^ 

may be rewritten in the form v'^. < Ur- Also 

Spec(M) n iSr - FniSr), Sr + Fn(Sr)) 7^ 

may be rewritten as 

Spec(M)n(z/,_i,z/;) ^0. 
Using (H) we now deduce that m, < z/'.. The proof that /i^ < /i^ < ''^r is similar. 

Improving the above result itcratively depends upon setting up the inductive hy- 
pothesis properly. Our hypothesis below allows the possibility of increasing the 
value of n from one stage of the iteration to the next, in other words of increasing 
the size of the test function space as the iteration proceeds. In practice one would 
choose a fairly small initial value of n and apply the lemma repeatedly until the 
iteration stops yielding improvements, and then increase n significantly and repeat 
the iteration with that value of n until improvements cease. 
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Theorem 8 Let /ir, /^^j ^r, J^'r, Sr and tr satisfy 
and 

Fn(tr) ^ tr 
l^r — Fji(tr) 

for 1 < r < R. Put z/q = (Tq and /i'r+i — (Ji?+i- Given (S) and (H), define /i'^, i/", 
s'j. and t'j. for 1 < r < R by 

Fn{tr) ~ A*r+1 ~ 

jjij. — t^ Fji(tj.^ 



Then 



for 1 < r < R. 

Proof Let 1 < r < R. The function 

(j){s) = Fn{s) - S + 

satisfies < for all s e R and 
Also 

(f){Sr) = Fn{Sr) - Sr + l^'r-l < K-1 " ^r-l < 0. 

Hence < s'^. < Sr- We next observe that 

- FniSr) <Sr-s', 
K = i^n(4) +S'r< F^{Sr) + Sr < V^. 

Spec(M) n (4 - 4 + 7^ 



implies 
Also 
is equivalent to 



Spcc(M) n {y'r-x.v^) ^ 
and implies that rur < I'r- The proof that fi'^ < /i" < rrir is similar. 
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At each stage in the implementation of the iterative process one has further consis- 
tency checks, namely the numerical procedure need not yield the inequality /i^' < v'^ 
for all r. If this does not happen then (H) must be false; that is the operator M 
must have more spectrum than assumed. 

Plotting the graph of a typical function Fn suggests certain results about the 
location of Sr and U] these are proved below. 

Lemma 9 Assume (S) and ([^. Given A G R, let X' denote the smallest point 
in Spec(M) greater than A. Assume A' > A + 4e and let s be the solution of 
F„(s) = s — A. Then 

A + A' A + A' ^ 

< s < h 2e. 

2 2 

Proof We combine the fact that the function (j){s) = Fn{s) — s + X is strictly 
monotone decreasing and continuous with the inequalities 

0((A + A')/2) > F((A + A')/2) - (A + A')/2 + A = 

and 

0((A + A')/2 + 2e) < F((A + A')/2 + 2e) + 2e-{X + X')/2 -2e + X = -2e<Q. 

Theorem 10 Suppose that M has a single eigenvalue m in {a, (3). Define s„ and 
tn by 

-fn('5n) Sfi tt 

Un = Sn + Fn{Sn) (= 2s„ - Ct) 
Fnij^n) /3 tn 

f^n tn Fnitn) ( 2tn 

Then 

fin < m < Un 

for all n > 1. Moreover 

lim s„ = (a + m)/2, lim t„ = (m + /?)/2 (14) 

and 

lim Un = lini Un = fn- 

n— ►oo n— >oo 

Proof The first statement follows from Lemma [7| with R = 1, Uq = a and H2 = 
The other statements are then consequences of Lemma IHl 
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6 The Zimmermann-Mertins Paper 



We have repeatedly referred to solving the equations = s—v and = n—t 

for various /x, v. This may be carried out by bisection or by Newton's method. By 
examining the paper of Zimmermann and Mertins, [201, one discovers a method 
which is probably better in most circumstances; see ()21|) below. 

We assume that (z/, /i) is an interval which contains a single eigenvalue m of M, 
but no other point of Spec (M). E.g. in the situation described in hypothesis (H), 
we have in mind tho choose v = Ur-i, = fJ^+i, ^ = f^ri for any r G {1, . . . , R}. 
Suppose for simplicity that z/ and /i are not eigenvalues of M. 



Let C denote an n-dimensional trial function space. To ensure that the 
exact eigenfunction corresponding to m is 'not too far' from we make the 
additional assumption 

{Mv,v) . {Mv,v) 
ly < max — — , mm — — < u. (15) 

(Without (fT3j) . the Lehmann method gives nothing, and the F„-method gives a 
lower bound for m less than z/, resp. an upper bound larger than n, which makes 
it more or less useless.) 



The right-definite Lehmann method may be applied to the situation in which only 
one eigenvalue is to be bounded. We follow fK)\ Theorem 4.2], where the spectral 
parameter p is once chosen to be /x, and once u. Let {vi, . . . , Vn} denote a basis of 
£n, and let r~ be the smallest eigenvalue of the n x n matrix eigenvalue problem 

{{Mvi - fiVi, Vk))i,k X = T{{Mvi - pvi, Mvk - pvk))i,k x (16) 

and r+ the largest eigenvalue of 

{{Mvi - uvi, Vk))i,k X = T{{Mvi - uvi, Mvk - i'Vk))i,k x. (17) 

Then 

uH <m<z/ + — . (18) 

r r+ 

(Note that < < r+, according to (US))). 

The following theorem states that the F„-method gives exactly the right-definite 
Temple-Lehmann-Maehly bounds, in the Zimmermann-Mertins setting. See also 
[71 Theorem 6]. 

Theorem 11 Let and as above, and let s and t denote the solutions of (cf. 
Lemma^: 

F^{s) = s-u, Fn{t)=p-t. (19) 
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Then, 



Hence 



s + = 1^ + ^, t- F^{t) =// + —. (20) 

. = . + ^, t = ^^ + ^. (21) 

Proof For p G R, define 

^2(P) := {{MVi - pVi, MVk - pVk))i,k=l,...,n, 

Ai{p) := {{Mvi - pvi,vk))i,k=i,...,n, Aq := {{vi,Vk))^,k=l,...,n■ 
StYa.ight^OTwa,Td calculations and use of (fT!I|l give, for all x G R", 

x^A2{s)x - Fn{sfx^Aox = x^A2{u)x + 2{u - s)x'^Ai{u)x 

+ [(^ _ ,)2 _ F^{sf]x^AoX 
= x^A2{u)x -2Fn{s)x^Ai{u)x. (22) 

By the definition of Fn{s), the left-hand side of ()22p is non-negative for all x, and 
zero for some specific x. This is therefore true also for the right-hand side, i.e., 

x'^Ai{l')x ^1 r u T^n\ rr^^ 

r . < TTTTTT for all X G R" \ {0}, 

a;M2(z/)a; 2F„(s) 

with equality for some specific x. Thus, l/(2F„(s)) is the maximum of the Rayleigh 
quotient of problem |T7jl. i.e., 

^ r+. (23) 



2F„(s) 

Using (fT^ again, we obtain 

S + Fn{s) = U + 2Fn{s) = 1^ + ^ ■ 

A similar method yields 

x'^A2{t)x - Fnit^x'^Aox = x'^A2{p)x + 2Fn{t)x^ Ai{p)x, 
and the same arguments as before yield 

^44^> fVt for all X G R" \ {0}, 

x^A2{p)x - 2F„(t) ^ ^ ^' 

with equality for some specific x, which implies 

1 



(24) 



2Fn(t) 
and thus 

t-F„(t) =/i-2F„(t) =/i + — . 

T 

The final statement of the theorem is obtaining by combining (fTIH) . and ()24|1 . 
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7 An Example in One Dimension 



Suppose that the operator H has a complete orthonormal set of eigenfunctions 
{0r}5^=i, whose corresponding eigenvalues Xr are all positive. Suppose also that 
the eigenvalues are listed in increasing order and tend to infinity. If we put 

iPr = K^'^D<Pr (25) 

then {il)r}'^=i forms an orthonormal basis in JV{. and the matrix elements of M with 
respect to this basis are 



The truncation method of computing the eigenvalues consists of choosing n large 
enough and then finding the eigenvalues of the matrix M„ obtained from M by 
considering only those coefficients for which r,s < n. This is equivalent to find- 
ing the spectrum of Mn = PnMPn, where P„ is the orthogonal projection onto 
Cn = hn{'?/'i, ■i/'n}. It is well known that every point of Spec(M) is the limit 
of eigenvalues of M„ as n — >■ oo. Spectral pollution arises from the fact that the 
converse is not true. 

As an example we consider the operator H acting in L^(0,7r) according to the 
formula 

dV 

and subject to Dirichlet boundary conditions. We then have 

4'nix) = \J 2 / sm{nx) 
and A„ = for n = 1,2, .... We take K to be the differential operator 

where b{-) is any non-negative bounded function on (0, vr). If D is the differentiation 
operator with domain iyQ^'^(0, vr) then H = D*D and K = D*BD, where B is the 
operator of multiplication by b{-) acting in L^(0, vr). The range of D is not L^(0, tt) 
but the subspace A4 of all functions / such that 

r fix) dx = 0. 
Jo 

This corresponds to the fact that in the present context ()25|) translates to 



i^n{x) = \J 2 / TV cos{nx) . 
Since n = 1, 2, ... these do not form a complete orthonormal set in L^(0, vr). 
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The bounded operator M is given for / G by 

M/ = 6/-(6/,^o)^o 

where iIjq{x) = vr"^/^. It follows by trace class perturbation techniques that M and 
B have the same essential spectrum. This equals the closure of the set of essential 
values of 6, or the closure of its actual set of values if it is piecewise continuous. 
However M may also have a further eigenvalue in any gap of its essential spectrum. 



Lemma 12 The isolated eigenvalues of M are of multiplicity 1 and are the nonzero 
solutions of the equation 

JO b[x) — A 

In particular if h is the characteristic function of the measurable set E C (0, vr) 
then M has only one isolated eigenvalue, namely \ = {-k — \E\) /t[ . 

Proof The task is to find all / G lv^(0, vr) which satisfy (/, ipo) = and 

&/-(&/,^o)V'o = A/ (26) 

where A is not in the essential range of b. We may normalize / by {bf, i/jq) = 1 
since if the LHS equals 0, we see from (I26|) that A lies in the essential range of b. 
We then insert 

/ = (6-A)^Vo 

into the normalization equation. The function S is analytic and has a positive 
derivative in each gap in the spectrum of M, so there is at most one eigenvalue in 
each gap. 

The entries of the matrix An = PnMPn of are given for 1 < r, s < n by 

2 

An,r,s = — i>{x) cos(rx) cos(sx) dx. 



7r JO 

The matrix entries of Dn = PnM*MPn are calculated using the formula 



{M*Mf,g)= / b{xyf{x)g{x)dx - {bf,^o) i^oM), 



vahd for all f, g E M. 

We consider the exactly soluble example 



bix) 



1 if < X < a 
otherwise. 



where < a < tt is arbitrary. The operator M then has two-point essential 
spectrum {0, 1} together with an eigenvalue m = {it — a) /ir of multiplicity 1. If 
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one puts a = 11/2 then m = 1/2; numerical calculations show that as n increases 
an increasing number of eigenvalues moving inwards from and 1 and 
apparently converging slowly to 0.5. The following table lists the four eigenvalues 
of An closest to 0.5. 



n 


Ai 


A2 


A3 


A4 


50 


0.0526 


0.4632 


0.5368 


0.9474 


100 


0.0691 


0.4750 


0.5250 


0.9309 


200 


0.0856 


0.4830 


0.5170 


0.9144 


400 


0.1019 


0.4884 


0.5116 


0.8981 



The use of the techniques of this paper not only avoids the problem of spectral 
pollution, but also provides accurate results for much smaller values of n. We start 
by putting a = tt/2 and n = 10. The function -Fio(A) has three local minima. 
They are given in the following table, in which the numbers are all accurate to 
eight decimal places. 

A Fio(A) 
0.00 0.00000679 
0.50 0.10049280 
1.00 0.00000679 

The central minimum is exactly equal to 1/2 as a result of the unusual symmetry 
for a = Tr/2. The hypothesis (H) is therefore valid with uq = 0, /ii = 0.39950720, 
Ui = 0.60049280 and /i2 = 1. An application of Lemma [7| now yields Si = 0.25 and 
ti = 0.75 to six decimal places, which confirms that the eigenvalue equals 0.5 to 
the same accuracy. 

We repeat the calculation for a = 1.0, in order to see whether the unexpected 
accuracy of the results are caused by the symmetry. Two representative graphs 
are shown in Figure 1. The three local minima of -F8(A) are now as follows, accurate 
to six decimal places. 

A Fs{X) 


0.697669 0.106699 
0.999313 0.026051 

The hypothesis (H) is therefore vahd with z/q = 0, /xi = 0.590970, z/i = 0.804368 and 
/i2 = 1. An application of Lemma [7| now yields Si = 0.340845 and ti = 0.840844 
to six decimal places. This implies that 

0.681688 < /i < 0.681690. 

The accuracy of the upper and lower bounds on the eigenvalue /i depends on the 
rate of convergence in ()14|) . which is closely associated to the rate of convergence 
in 

lim = lim F„(t) = F{t) 
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where s = m/2 and t = (m + l)/2. For reasons which we do not understand 
the rate of convergence at these two points is much faster than at other points in 
(0, 1). This phenomenon is illustrated in the following table, in which u and v are 
'general' points. 

Putting a = 1.0 we have m = 0.681690, s = 0.340845 and t = 0.840845. We also 
put u = 0.5 and v = 0.75. 



n 










2 


0.340873 


0.269294 


0.254218 


0.433012 


4 


0.340845 


0.162575 


0.238589 


0.146875 


6 


0.340845 


0.159219 


0.221431 


0.137537 


8 


0.340845 


0.159156 


0.210113 


0.116463 


10 


0.340845 


0.159155 


0.207660 


0.108933 


100 


0.340845 


0.159155 


0.187173 


0.077840 


oo 


0.340845 


0.159155 


0.181690 


0.068310 



The final row lists the values of F{-). 

The surprisingly rapid convergence of Fn{s) and Fn{t) does not occur only for the 
two values of a presented. A further computation indicates that if we put s = m/2 
then 



for all < a < 2.3. For larger a the difference Fiq{s) — F{s) increases steadily. 
We also have 



s = F{s) < F5o(s) < F{s) + 10-^ 

for all < a < 2.9. 

8 An Example with Variable Coefficients 



It might be thought that the cause of the spectral pollution in the above examples 
was the fact that the basis chosen was badly adapted to the operator. 

In this section we consider the same type of example as in the last, but using 
a different type of basis, constructed using the Krylov subspace method. This 
example does not display spectral pollution, but this is not because the associated 
matrix is tridiagonal. A recent example of Denisov and Simon shows that spectral 
pollution can occur in the following sense even for tridiagonal (Jacobi) matrices, 
jlOj . They construct a bounded, self-adjoint, tridiagonal operator A acting on 
/^(Z"*") with spectrum [—5, —1] U [1,5]. The existence of a gap in its spectrum is 
not a priori obvious. For certain values of n, the truncation An of A to the space 
Cn of sequences with support in {1,?t,} has a single eigenvalue in (—1, 1). The set 
of limit points of these eigenvalues is equal to the whole interval [—1, 1]. 

We replace the interval (0, vr) by (—2, 2) and put 



s = F{s) < Fio(s) < F{s) + 10 
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The essential spectrum of b, and of M, is thus [a" — 2P~,a~] U + 2/?"'"]. 

Assuming that > and a~ there is an eigenvalue A in the spectral gap 

This is evaluated by solving 

as in Lemma IT^ 

We next write down the basis which we will use. Let {Vn}'^=i denote the Legendre 
polynomials, and define 

SO that 

for all n > 0, where 

n 

In 



Then define on (-2, 2) by 



Finally define ipn for all n E Z by 



Qn{x - 1) if ,T > 

if X < 0, 

if X > 

Qn{x+1) ifx<0. 




if n > 1 

= { Q\n\ if ^ < -1 

iQo+Qo)/V2 ifn = 0. 
It is easy to check that {ipn}nez is a complete orthonormal set in 

M = {fe 2) : £ fix) dx = 0}. 

Let Pn be the orthogonal projection of onto = \m{ip_n, ■■■,4'n}- Then it 
follows from standard properties of the Legendre polynomials that M(£„) C Cn+i- 
Hence M has a tridiagonal matrix. Routine calculations yield 

a+ + P+ if r > 1 

Mr,r = { a~ - if r < -1 

{a+ + (3+ + a- - f3-)/2 ifr = 0. 



Also 



{ /3+7r+i ifr > 1 

/3+7i/y2 if r = 

/3-7i/y2 ifr = -l 

[ r7M ifr<-2. 
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All other coefficients of M vanish. We observe that Mr^r and M^^r+i have different 
limits as r — » ±00. Operators whose coefficients have different asymptotics on 
the left and right have an interesting spectral and scattering theory, studied, for 
example, in 0. 

We finally compute Bn{X) from (fTT|) . obtaining the formula 

5„(A) = (A„-A/)2 + E„ 
where En = PnM{I - P„)MP„. Explicitly we have 

( {p-'jn+iY if r = s = -ra 

En,r,s = I [P^ln+lf H r = S = U 

[ otherwise. 

We carried out numerical calculations for the case = 1, a~ = —1, f3~^ = 40, 
/3~ = 30, n = 50, so that dim(£„) = 2n + 1. The essential spectrum of M is 
[—61, — 1] U [1, 81] and by solving pTjl we find that M also has a single eigenvalue 
/i ~ 0.453261434. Plotting the function F{\) on [0, 1] shows that it has a well- 
defined local minimum near this point and nowhere else in the interval. From the 
values 

F5o(0.23) = 0.22326196 
F5o(0.72) = 0.26673868 

we obtain the enclosures 

0.72 - F5o(0.72) <fx< 0.23 + F5o(0.23) 

or, more explicitly /i = 0.45326^°. Of course all these computations could be 
carried out using interval arithmetic, but we did not do so. These upper and lower 
bounds on A compare with the completely uncontrolled values for A obtained by 
the truncation method. 

n fjL 
10 0.50275311 
20 0.45504606 
50 0.45326153 
75 0.45326143 
100 0.45326143 

In this example we have found that by evaluating -F5o(A) at only two points, we 
obtained upper and lower bounds on fi which were of comparable accuracy to the 
uncontrolled value obtained by truncation. The full plot of Fn{X) gives valuable 
geometrical insights, but may be carried out for a smaller value of n than is used 
for the final bounds; this is, of course, much faster. 

Acknowledgments We would like to thank M Levitin and E Shargorodsky for 
valuable discussions. 
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Figure 1: Fs{\) and Fioo(A) for a = 1.0 
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